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The complex scaling method, which consists in continuing spatial coordinates into the complex 
plane, is a well-established method that allows to compute resonant eigenfunctions of the time- 
independent Schrodinger operator. Whenever it is desirable to apply the complex scaling to inves- 
tigate resonances in physical systems defined on numerical discrete grids, the most direct approach 
relies on the application of a similarity transformation to the original, unsealed Hamiltonian. We 
show that such an approach can be conveniently implemented in the Daubechies wavelet basis set, 
featuring a very promising level of generality, high accuracy, and no need for artificial convergence 
parameters. Complex scaling of three dimensional numerical potentials can be efficiently and ac- 
curately performed. By carrying out an illustrative resonant state computation in the case of a 
one-dimensional model potential, we then show that our wavelet-based approach may disclose new 
exciting opportunities in the field of computational non-Hermitian quantum mechanics. 



I. INTRODUCTION 

Resonant states, sometimes referred to as metastable 
quasi-bound states, are eigenstates of the time- 
independent Schrodinger operator (atomic units are used 
throughout this work), 



H[r\ = f[r\ + V[r\ 



-\v 2 + V{r), f=(x,y,z). (1) 



These solutions exhibit complex- valued eigenvalues, the 
(negative) imaginary part being inversely proportional to 
the resonance finite lifetime. Several equivalent defini- 
tions of resonant states can be formulated, each one sug- 
gesting a different computational approach [1]. For in- 
stance, resonant states can be either regarded as poles of 
the scattering matrix, or as solutions of the Schrodinger 
equation under the condition that only outgoing asymp- 
totes exist (Siegert boundary conditions [2]). 

When present, resonant states are not to be found 
within the Hcrmitian spectrum of H, nor the correspond- 
ing eigenfunctions are elements of the Hilbert space of 
square-integrable functions. The fundamental works of 
Aguilar, Balslev, Combes and Simon [3-5], proved that 
resonant states show up as discrete states in the spectrum 
of a non-Hermitian Hamiltonian operator obtained by 
scaling spatial coordinates by a complex factor, namely: 



-> H[fe i8 } = H e 



(2) 



As soon as the complex scaling angle 9 is greater than 
the critical angle 9 C = axcta,n[\Ei/(E R — E th )\]/2, E t h = 
V(oo) being the threshold energy, the resonant state cor- 
responding to the eigenvalue E = En — iEj, Ej > 0, 
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appears as a discrete eigenstate of Hg , with regular eigen- 
functions tpg(r) — > as f — > oo (see Figs. 3 and 4 for an 
illustration) . 

A conspicuous number of articles contributed to defin- 
ing a non-Hermitian formalism in quantum mechanics. 
This formalism generalises all the fundamental ingredi- 
ents of traditional quantum mechanics (Rayleigh-Ritz 
variational principle, Hellmann-Fcynman theorem, ...), 
rephrased in terms of a modified inner product ("c- 
product" or bi-orthogonal product) defined as follows: 



if\g) = {r\g) 



(3) 



Here, (f\g) represents the usual inner product in the L2 
Hilbert space. For more comprehensive information, wc 
refer to e.g. Refs. 6-12. 

Obtaining the complex-scaled version of the kinetic op- 
erator T is straightforward, since f g = e- i2S f. However, 
the evaluation of the complex rotated potential operator 
Vg is less obvious, especially if the potential is not known 
in closed form. For instance, this is always the case within 
self-consistent numerical simulations. To this end, differ- 
ent options have been proposed in the literature. Sum- 
marizing, they can be sketched out as follows: either the 
numerical potential is fitted with polynomials up to some 
order and then complex scaling is applied to the matrix 
elements in some finite basis set representation [13], or 
the unsealed potential is represented on a complex-scaled 
basis set [14, 15]. 

Complex-scaled Hamiltonians can also be obtained 
through a similarity transformation [16], without having 
to carry out the explicit evaluation of the original Hamil- 
tonian in a complex configuration space (see e.g. Rcf. 
6, p. 153). Within this approach, the complex-rotated 
Hamiltonian is obtained as 



Hg — SgHSn 



(4) 



where the linear complex scaling operator Sg is such that 



2 



fg = e-' l29 f and 

[SgV](r) = V(re w ) . (5) 

The mapping between the c-orthonormal eigenfunc- 
tions of H and Hg reads = e'^^Sg 1 \ip g ). In the 
case of potentials for which a Taylor expansion is well- 
defined, the similarity transformation is strictly equiva- 
lent to the analytic continuation of the potential into the 
complex plane. 

Different options for the operator Sg are possible, as 
extensively discussed in Ref. 16. We consider here the 
simplest formulation, namely: 

Sg = e eW , with W = i(r-f )- V, (6) 

r*o = (a^yo^o) being the fixed-point of the rotation. 
Clearly, Sg 1 = S-g due to the group structure of the 
complex scaling transformation. 

Although conceptually simple, implementing the com- 
plex scaling method through a similarity transformation 
in a discrete basis set is not a trivial task. A remarkable 
attempt can be found in Ref. 17, where the similarity 
transformation operator was represented in the sine dis- 
crete variable representation (DVR) basis set. The accu- 
racy of the method was ultimately proven to be somewhat 
sensitive to an artificial convergence parameter, depend- 
ing on the value of 9 and on the basis set size. 

In this study, we show that Daubechies wavelets enjoy 
a number of features that allow obtaining a reliable rep- 
resentation of all the operators invoked by the complex 
scaling method, with no need of artificial convergence 
parameters. 

The plan of the paper is as follows: after briefly re- 
calling the wavelet formalism, we show in Sees. II and 
III how to represent the Schrodinger Hamiltonian and 
the complex scaling operator in a wavelet basis set. In 
Sec. IV, we review some of the standard results of non- 
Hermitian quantum mechanics, which involve the com- 
plex virial theorem and eventually allow to assess the 
rate of convergence towards exact bound and resonant 
states. In Sec. V, we present the results that we ob- 
tained when applying our wavelet-based algorithm to a 
number of model potentials, including a multi-centered 
one-dimensional (ID) potential and a three-dimensional 
(3D) potential. We will show that the algorithm exhibits 
excellent accuracy and numerical stability. As a further 
confirmation, we will extensively discuss the computation 
of the lowest-lying resonances exhibited by a ID model 
potential, so as to prove that the complex virial theorem 
turns out to be fulfilled down to machine precision. Our 
conclusions are summarized in Sec. VI, along with an 
outlook on future developments. 

II. OPERATORS IN WAVELET BASIS SETS 

The adoption of wavelets [18] as a basis set in the 
present context is advisable for a number of reasons. 
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FIGURE 1: Least asymmetric Daubechies wavelet family of 
order 2m = 16. Note that both the scaling function <fi(x) 
and the wavelet if>(x) are different from zero only within the 
interval [1 — m,m]. 

First, Daubechies wavelets [19] present several properties 
which make them suitable for the numerical simulation of 
isolated systems: they form a systematic, orthogonal and 
smooth basis, which is localized both in real and Fourier 
spaces and allows for adaptivity. Second, we will show 
that such functions will meet both the requirements of 
precision and localization found in many applications of 
the complex scaling method. While referring the reader 
to Ref. [18] for an exhaustive presentation of how wavelet 
basis sets can be used for numerical simulations, we here 
summarize the main properties of Daubechies wavelets, 
with a special focus on the representation of the objects 
(wavefunctions and operators) involved in the present ap- 
proach. 

Every wavelet family comprises a scaling function <f>, 
and a second function ip properly called wavelet. Fig. 1 
illustrates the least asymmetric Daubechies wavelet fam- 
ily of order 2m = 16, the basis set which is used in the 
present study. These functions feature a compact sup- 
port [1 — to, to] and are smooth, therefore localized in 
Fourier space as well. 

A basis set is simply generated by the integer translates 
of the scaling and wavelet functions, with arguments mea- 
sured in units of the grid spacing t. For instance, a ID 
domain of extension L, centered at x = 0, can be spanned 
by the following set of N scaling functions, 

(xfa) = &(z) =-U(f"i),i= 'N/2, . . . , N/2 , 

(7) 

where £ = L/(N — 1) is the (uniform) grid spacing. The 
basis set can be completed by the addition of the trans- 
lates of the wavelet functions ipi. These functions form a 
orthogonal basis set: 

(&|&>=*« = MI^>, W>^)=0. (8) 

In three dimensions, a wavelet basis set can easily be 
obtained as the tensor product of one-dimensional basis 
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functions, combining wavelets and scaling functions along 
each coordinate of the Cartesian grid (see e.g. Ref. 20). 

The most important feature of any wavelet basis set is 
related to the concept of multiresolution. Such a feature 
builds upon the following scaling eguations (or "refine- 
ment relations"): 

</)(*) = V2]TM(2z-j); (9) 

3 

i>{x) = ^^gjtjtQx-j), (10) 

i 

which relate the wavelet representation at some resolu- 
tion to that at twice the given resolution, and so on. Ac- 
cording to the standard nomenclature, the sets of the hj 
and gj = (— l) J /ii_j coefficients are called low- and high- 
pass filters, respectively. A wavelet family is therefore 
completely defined by its low-pass filter. In the case of 
Daubcchies-2m wavelets, j G [1 — to, to]. 

The representation f(x) of a function in the above de- 
fined basis set is given by: 

N/2 N/2 

f( x ) = CiM x )+ ^2 diipi{x), (11) 

i=-N/2 i=-N/2 

where the expansion coefficients are formally given by 
Ci = (<f>i\f), di = (VM/)- Using the refinement equations 
(9) and (10), one can map the basis appearing in Eq. (11) 
to an equivalent one including only scaling functions on 
a finer grid of spacing £/2. 

The multiresolution property plays a fundamental role 
also for the wavelet representation of differential opera- 
tors. For example, it can be shown that the exact matrix 
elements of the kinetic operator, 

Ti-j = -\\ Ax^{x)d 2 ^{x) , (12) 

equal the entries of an eigenvector of a matrix which 
solely depends on the low-pass filter (see e.g. Ref. 18). 
The construction of the kinetic operator T in a wavelet 
basis set does not require any sort of approximation or 
numerical integration. In the next section, wc will show 
that the matrix elements of W can be calculated in a 
similar way. 

The discretization error due to Daubechies-2m 
wavelets is controlled by the grid spacing. Daubcchics- 
2m wavelets exhibit to vanishing moments, thus any poly- 
nomial of degree less than to can be represented exactly 
by an expansion over the sole scaling functions of order 
to. For higher order polynomials the error is 0(£ m ), 
i.e. vanishingly small as soon as the grid is sufficiently 
fine. Hence, the difference between the representation 
of Eq. (11) and the exact function / decreases as £ m . 
Among all the wavelet families, Daubechies wavelets fea- 
ture the minimum support length for a given number of 
vanishing moments. 



Given a potential V known numerically on the points 
{xk} of a uniform grid, it is possible to identify an ef- 
fective approximation for the potential matrix elements 
Vij = (<f)j\ V \<f>i). It has been shown [20, 21] that a 
quadrature filter {w^} can be defined such that the ma- 
trix elements given by 

Vn = {4>j\ V |&) = w h -i V(x k ) uJk-j , (13) 

k 

yield excellent accuracy with the optimal convergence 
rate 0(£ 2m ) for the potential energy. The same quadra- 
ture filter can be used to express the grid point values of 
a (wave) function given its expansion coefficients in terms 
of scaling functions: 

f{x k ) = Y d CiUk-i + G{i m ); (14) 

i 

c i = J2f(x k )LJ k . i + 0(£ m ). (15) 

k 

As a result, the potential energy can equivalcntly be com- 
puted either in real space or in the wavelet space, i.e. 
(f\V\f) = E k f(xk)V(x k )f(x k ) = EijdVijCj. The 
quadrature filter elements can therefore be considered 
as the most reliable transformation between grid point 
values f(x k ) and scaling function coefficients Cj, as they 
provide exact results for polynomials of order up to m — 1 
and do not alter the convergence properties of the basis 
set discretization. The filter {ujk} is of length 2to and 
is defined unambiguously by the moments of the scaling 
functions (which in turn depend only on the low-pass 
filter) [22, 23]. 

Using the above formulae, the Hamiltonian matrix 
Hij — Tij + Vij can be constructed. Note that, in con- 
trast to other discretization schemes (finite differences, 
DVR, plane waves, etc.), in the wavelet basis set neither 
the potential nor the kinetic terms have diagonal repre- 
sentations. Instead, H is represented by a band matrix 
of width 4to — 1. 



III. COMPLEX SCALING VIA SIMILARITY 
TRANSFORMATION 

The complex scaling operator Sg being separable along 
the three spatial directions, and the wavelet transform of 
3D objects being given by the tensor product of three 
ID wavelet transforms, it is sufficient to build the com- 
plex scaling operator suitable to working on ID input 
functions. 3D complex-scaled potentials can thus be ob- 
tained through the threefold application of the same ID 
operator, each time acting along a different direction: 

V(xe w ,ye w lZ e ie ) = [§$§»§$V](x,y,z). (16) 

In what follows, we drop the superscript denoting the 
spatial direction and focus on the construction of the ID 
complex scaling operator, taking 

W = i(x-x )^-. (17) 
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FIGURE 2: Spectrum of the operator W - cf. Eq. (17) - for 
different sizes N of the basis set. Knowing that the spectrum 
of W is symmetric around the imaginary axis, only the right 
half-plane is shown. 



Owing to the scaling properties of wavelets and scaling 
functions, the matrix elements Wij of the operator W 
can be computed exactly. The final result reads as 



W, 



Wi-4 + 



i+j 



x 



di — q 



(18) 



The definitions of the filter coefficients {wj} and {dj} are 
deferred to Appendix A, along with the details of their 
computation. Let us note that, within the Daubechies- 
2m family, Wi j = for any pair {i, j} such that \i — j\ > 
2m — 2. Moreover, being W scale- invariant, its spectrum 
depends only on the size N of the basis set. Fig. 2 shows 
that the eigenvalues of W span a range which reaches a 
maximum that is < N. 

Once the exact wavelet representation of the operator 
W is available, the operator Sg can be represented by its 
spectral decomposition over the eigenvalues and eigenvec- 
tors of W: 



N 



So 



jnv 



= E e * Ai i^xv 



Ail 



(19) 



real and reciprocal space. Actually, since the exact eigen- 
function of W corresponding to the (complex) eigenvalue 
A is given byi/j \ (x) cx (x— x )~ lA , the real part of A can be 
deemed as a wavenumber in a logarithmic configuration 
space. Given that xd x = —(1 + pd p ), the same behav- 
ior would hold also in reciprocal space, ij)\(p) oc p lA_1 . 
Smooth and localized functions are therefore likely to ex- 
hibit non-trivial projections only onto the lowest-lying 
cigenstates of W. 

At the same time, for two functions f{x) and g{x) = 
f(x) + c that differ by a global constant c € C, the equal- 
ity Sgg = Sgf + c is expected to hold for any 9. The nu- 
merical representation of the operator Sg should namely 
behave as the identity on a generic constant vector |c). 
A numerically stable implementation should therefore be 
such that (ipXilV) = (ipx t \V + c) for any Aj ^ 0. Such a 
requirement can be enforced by imposing suitable bound- 
ary conditions on the matrix elements of W , as demon- 
strated by the following result, which involves the ?i-th 
moment of a generic left-cigenfunction ip\( x ) °f W (as- 
suming x = 0, Xipx(x) ee ipx(xjw = id x [xipx(x)]): 



L/2 

(A + in) / dxip\(x)x n = 

-L/2 

I ^y +1 [ML/2) + (-l) n M-L/2)} . (20) 



In particular, for n — the boundary term on the right 
hand side becomes zero if ip\(L/2) = —ijj\(—L/2), and 
(ipx\c) = for any A ^ 0. In other words, the only left- 
eigenvector of W that has a non-trivial projection on a 
constant, is the one corresponding to the zero eigenvalue. 
Such boundary conditions can be imposed by setting 



Wi. 



j±N 



Wi-4 + 



i+j 



X 



(21) 



The zero eigenvalue can actually be found in the spec- 
trum of the W operator that is so obtained (see Fig. 2), 
the pair matching right-eigenvector being identically con- 
stant. 



where A-s is the eigenvalue corresponding to the eigen- 
vector \ip\i)- However, special attention has to be paid 
when using Eq. (19) as is, since numerical instabilities 
are likely to occur depending on the maximum value 
attained by Re(#Ai) (in other words: on the condition 
number of the operator Sg). For instance, if the com- 
putation is carried out in double-precision floating point 
arithmetic (by far the most typical case), the relevant 
terms of the spectral decompositions (19) should lie in 
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g\i 



< 10 15 . In order for the com- 



IV. THE COMPLEX VIRIAL THEOREM 

As a result of the Balslev-Combes theorem [4], bound 
and resonant eigenstates are expected to be stationary 
with respect to variations of the complex scaling an- 
gle 9. Moreover, the complex analog of the Hcllmann- 
Feynmann theorem [6] allows to relate the variation of 
the i-th eigenvalue E^ 1 ' with respect to 9 (considered as a 
variational parameter) to the quantum expectation value 



plex rotation to be reliable, the scalar product (tpxAV) of dHg/d9 on the corresponding eigenstate \tp\ ): 
should be non-zero only for a subset {A^} for which |e' 



ft A, 



falls within the same range. This is however likely to hap- 
pen, provided the potential V is sufficiently localized in 
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FIGURE 3: Left: the results of the application of the wavelet-based complex scaling method to the potential (29) for different 
complex scaling angles. In the bottom row, we plot the absolute value of the difference between the exact complex-scaled 
potential V(xe l9 ) and that obtained through similarity transformation in the wavelet basis set, [SeV](a;). Right: the low-energy 
portion of the spectrum of the corresponding complex-scaled Hamiltonian. One can note that continuum states align along 
rays (dot-dashed straight lines in the plot) that are rotated by 20 with respect to the real axis, while bound and resonant states 

can be identified by the fact that they do not depend on 6. 



being 



E 



(23) 



Given the definition of Hg in Eq. (4), it is easy to show 
that 



dH g 
86 



= [W,H e ], 



(24) 



and to rephrase the requirement of ^-independence in 
terms of the complex virial theorem (CVT): 

— - {^WW^HeW^) =0. (25) 
As shown e.g. in Ref. 24 it is possible to prove that 

[W,H e ] = -2i(fg + U e ), (26) 

where 



2 dx 2 
Eventually, the virial theorem reads as follows: 



Tg + Ug 



When using finite basis sets, the condition (28) may not 
be satisfied, or be satisfied to a great extent only within 
certain intervals of the angle 9. Such an occurrence is ex- 
tensively reported in the literature [25]. Resonant eigen- 
values have been shown to move along the so-called 6- 
trajectories in the complex plane as 8 is varied, the fea- 
tures of which (cusps, in particular) can be used to assess 
the level of convergence towards stationary points, and 
might hint at the selection of the optimal complex scaling 
angle [13, 26, 27]. 

The evaluation of dE^ / d8 thus represents an assess- 
ment of the degree of convergence towards true bound 
and resonant states in numerical computations with fi- 
nite basis sets. A closer look at Eq. (27) reveals that no 
additional effort is required to obtain the operator U in 
the wavelet basis set, once W has been generated. More- 
over, the complex-scaled counterpart Ug can be obtained 
by acting on U with the complex scaling operator Sg, 
exactly as if U were a generic potential operator. 



(27) V. ILLUSTRATIVE EXAMPLES 

In order to validate our wavelet approach to complex 
scaling via similarity transformation, we first discuss the 

(28) 

application of our method to a specific ID model poten- 
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FIGURE 4: Complex-rotated c-normalized eigenfunctions corresponding to the five lowest-lying eigenstates (from left to right: 
one bound state and four resonant states), computed at 9 = {0,0.2,0.3,0.4} rad. Solid (dashed) lines correspond to the real 
(imaginary) part of the eigenfunctions. Dot-dashed lines represent the potential (rescaled by a numerical factor), which is 
superimposed so to highlight the localization of the complex-rotated eigenfunctions inside the interaction region. 
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FIGURE 5: Measurement of the fulfillment degree of the 
complex virial theorem, Eq. (28). The result of the evalua- 
tion of Eq. (28) is shown as a function of the grid spacing of 
the simulation domain, for the five lowest-lying states of the 
complex-scaled Hamiltonian at 9 — 0.4 rad obtained through 
similarity transformation in the Daubechies-16 wavelet basis 
set. The potential is that of Eq. (29). 



tial. We chose to work on the following potential, 
'1 o 



V(x) 



-x — a e 



-\x A 



bh 



(29) 



with a = 0.8, A = 0.1 and E t h — 0, because the same 
potential, with E t h = 0.8, has been thoroughly studied 
in the literature, hence we could cross-check our results 
against the already published ones. Let us point out that 
we do not aim to provide new information concerning the 
resonances exhibited by the model potential (29), instead 
to show that 1. the wavelet-based complex scaling trans- 
formation is extremely accurate; 2. the finite wavelet ba- 
sis set computation of resonant states yields very good 
figures of the fulfillment of the CVT. 

The potential (29), which models diatomic molecular 
pre-dissociation resonances, has already been studied in 
a number of publications, namely in Refs. 27, 28 (adopt- 
ing a Gaussian basis set), in Refs. 29, 30 (employing 
Weyl's theory and the direct numerical integration of 
the complex-rotated Schrodingcr equation under Siegert 
boundary conditions), in Refs. 31, 32 (employing a non- 
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FIGURE 6: The results obtained for Ho with the potential 
as from Eq. (29), with a focus on the low-energy end. The 
grid spacing was set equal to 0.075 (cf. Fig. 5). The offset 0.8 
(= —V(x = 0)) was added to the eigenenergies in order to 
avoid ambiguities in the ordering of states, even in terms of 
the absolute values of the associated eigenenergies. Vertical 
dashed lines are superimposed on the eigenenergies of the five 
lowest-lying discrete states. 



Hcrmitian generalization of the Milne's method), in Ref. 
33 (using a complex absorbing potential and plane waves 
as basis functions) and in Refs. 15, 34-37 (using the sine 
DVR). More specifically, in Refs. 36, 37 the potential is 
complex-scaled while the basis functions are not; in Refs. 
15, 34 and 35 the basis functions are scaled, but not the 
potential. The model potential (29) exhibits only one 
bound state (at E — E^ ~ —0.3), and several resonant 
states. 

The outcome of the complex scaling operation applied 
to the potential (29) is shown in Fig. 3, together with 
a comparison against the exact complex-scaled potential. 
An excellent agreement is found throughout the spatial 
range for all the values of 9. 

The full diagonalization of the complex-scaled Hamilto- 
nian obtained through similarity transformation provides 
a further confirmation of the reliability of the wavelet- 
based complex scaling operation. To this effect we show, 
in the right panel of Fig. 3, the low-energy region of the 
Hg spectrum and, in Fig. 4, the eigenfunctions of the 
five lowest-lying states (one bound state and four reso- 
nances) . It can be observed that resonant eigenfunctions 
are well localized inside the interaction region and well- 
behaved throughout the spatial range, pretty much like a 
bound state. The localization increases while increasing 
the complex scaling angle (the effect being more evident 
for the resonances at higher energy). Also shown are the 
eigenfunctions sorted out from the spectrum of Hg = o by 
looking up the eigenvalues which were the closest to the 
resonant eigenvalues computed at 8 > 0. 

In Fig. 5, the rates of convergence towards true reso- 
nant states (namely, the degree of fulfillment of the CVT) 
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FIGURE 7: The result of the wavelet-based complex scaling 
as applied through similarity transformation to an illustrative 
multi-centered potential V(x) = 0.25 exp [-0.1(2; - 7.5) 2 ] + 
0.5exp[-0.1(x + 5) 2 ] - exp(-0.1a; 2 ). Left: real (top) and 
imaginary (center) part of the complex-scaled potential at 
9 = 0.4 rad, and comparison (bottom) against the exact re- 
sult. Right: low-energy portion of the spectrum of the cor- 
responding non-Hermitian Hamiltonian. By virtue of the 6- 
independence one can identify a pair of bound states and a 
number of resonances. 



are displayed. The rates improve upon decreasing the 
grid spacing, as a consequence of the systematicity of 
Daubechies wavelets. In the case of a sufficiently fine 
grid, the figures of the violation of the CVT become as 
little as < 10~ n , for all five lowest-lying states. Such 
figures have to be compared with those referring to con- 
tinuum states (> 10 -2 , sec Fig. 6). The gap of several 
orders of magnitude between bound/resonant states and 
continuum states allows to easily distinguish the former 
from the latter within the Hg spectrum. It is also inter- 
esting to note how the measure of the violation of the 
CVT changes vs 9 for the different resonances. 

We wish now to discuss the results obtained upon ap- 
plying the method to a multi-centered one-dimensional 
potential (Fig. 7) and to a 3D potential (Fig. 8), with the 
twofold aim of providing evidence that the method also 
works in three dimension and that no special restrictions 
have to be imposed on the symmetry of the potential. 
While referring to the figure caption for further informa- 
tion on the shape of the two potentials, we wish to remark 
that the wavelet-based complex scaling method proves to 
be, once again, very accurate. 

We look forward to applying the wavelet-based com- 
plex scaling method to genuine ab initio potentials, as a 
first crucial step towards resonances in realistic 3D (iso- 
lated) systems. 
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FIGURE 8: The result of the wavelet-based complex 
scaling applied to a fictitious 3D potential, V{r) — 

ELo afce ~ /3fc|r_rfc|2 with a ° = ~ 3 ' ai = 2 ' a2 = 4 ' ffl3 = 3 > 
fa = 0.1, /3i = /3 2 = = 0.2, r = 0, n = 3(0,-1,1/2), 

r 2 = 3/2(— v/3, 1,0), r 3 = 3/2(1, V^,0) and 9 = 0.3 rad: 2D 
surfaces at z = representing the real (top), imaginary (cen- 
ter) part of [SeV](r) and the discrepancy (bottom) from the 
exact result. Note that, for better visibility, the opposite of 
the imaginary part is plotted. 



VI. CONCLUSION 

We have shown that the adoption of wavelets as a ba- 
sis set allows to implement the complex scaling method 
via similarity transformation in a rigorous and effective 
way. Since wavelets display well-defined transformation 
properties upon rescaling of the spatial coordinates, they 
are especially well-suited for the discretization of scale- 



invariant operators such as the complex scaling generator 
W . In addition, the localization of wavelets (both in real 
and reciprocal space) allows to represent complex scaling 
transformation through the spectral decomposition over 
the eigenstates of a matrix that is inherently band-like. 
As a consequence, no artificial convergence parameter has 
to be introduced. 

By testing the method on a host of ID and 3D model 
potentials, we were able to prove that our implementation 
is also very accurate. This general technique, starting 
from the numerical values {f(xk)} of a generic function 
f(x) on a (uniform) real space grid, computes the values 
{f(xke lS )} by performing a similarity transformation in 
the Daubcchies wavelets basis set. The output may be 
further processed in any other basis set of choice, e.g. for 
the extraction of resonant states of the complex-scaled 
Hamiltonian. 

By carrying out a resonant state computation in the 
case of a one-dimensional model system, we were able 
to demonstrate that our approach exhibit excellent con- 
vergence rate in terms of the size of the basis set. The 
violation of the complex virial theorem can become com- 
parable to the machine precision, hence totally negligible. 

Resonant states are supposed to provide a concise de- 
scription of one-particle excitations, meaning that a cer- 
tain energy range ought to be spanned by a few reso- 
nant states rather than a virtually infinite set of contin- 
uum states. It has been claimed, for instance, that 11 one 
can construct the Green function simply in the form of 
a sum over the Siegert states, avoiding the annoying in- 
tegral over the continuum? (citation from Ref. 38; see 
also Refs. 1, 39). Several other investigations indicate 
the possibility of describing scattering cross sections, op- 
tical absorption spectra et cetera by means of a collec- 
tion of purely discrete (bound and resonant) states [40]. 
Nevertheless, the very numerical computation of resonant 
states is a formidable challenge. In our study, we showed 
that wavelets allow obtaining a complex-scaled Hamilto- 
nian which can contain resonant states in its spectrum (if 
any). Irrespectively of the basis set in which one wishes 
to further the computation, in realistic numerical sim- 
ulations the size of the Hamiltonian matrix can be so 
large that the full diagonalization becomes practically un- 
affordable. Moreover, one has to cope with the fact that, 
within the framework of non-Hermitian quantum mechan- 
ics, the variational principle states that "true" states are 
only stationary points along the variational trajectories 
(instead of minima). Hence, methods that are already 
in use for the computation of bound states cannot be 
deployed as such. 

A crucial achievement toward the numerical computa- 
tion of resonant states should regard the development 
of clever techniques for their direct extraction, allowing 
to focus on selected sub-regions of the spectrum of the 
non-Hermitian Hamiltonian (like, for instance, the filter 
diagonalization method [41]). Work is in progress in this 
direction. We also look forward to carrying out computa- 
tions in the frame of many-body perturbation theory with 
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the inclusion of resonant states [42]. hoping to be able to 
shed new light on the theoretical description of optical 
and electronic excitations. To this respect, the present 
study represents a first essential step, since it opens up 
the possibility of obtaining the complex-scaled version of 
Hamiltonians based on 3D ab initio potentials. 
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Appendix A: Representation of the W operator in a 
wavelet basis set 



Computing the wavelet representation of the operator 
W - cf. Eq. (17) - amounts to evaluating the following 
matrix elements: 



on </, and a non-translation-invariant part. Eventually, 



w m = 

r r rn.n — 


/ 6.x 4> m (x)W 4> n {x) ; 

'R 


(Ala) 


w< 01 > = 

r r m,n — 


/ dx m (x)Wip n {x) ; 

'R 


(Alb) 


= 

r r m,n — 


/ dxip m (x)W(f> n (x) ; 

'R 


(Ale) 


W {X1) = 


/ dxip m (x)Wil) n (x) , 


(Aid) 



where 4>i(x) = <fi(x — i) and ipi(x) = %p(x — i) are the 
scaling function and the wavelet belonging to a wavelet 
family, respectively. In our specific implementation, we 
adopted Daubechies wavelets, for the reasons outlined in 
Sec. II. Introducing the following variables, 



m + n 



J = 



(A2) 



and switching from x to x' = x—I, we can split Wm°n into 
a translation-invariant part, which indeed depends only 



J 



W m = i 
where 



(00) 



m + n _ _ \ (oo) 
Xq I a m _ n 



, (A3) 



w„ (00) = / dx4>(x-i/2)x^-(j){x + i/2) (A4) 

ox 



and 



,(oo) 



dx(f>(x — i)—(j>(x) — — 



(A5) 



The last equality is a consequence of the scaling func- 
tions and wavelets' compact support. The results for 
Wm}n, Wm°n, Wrn.n follow analogously, provided that 
in the preceding equations we change </>(cc) — > ^(x) and 
— >• 1 according to Eqs. (Al). 

The entries of the vector d^ 00 ^ can be computed follow- 
ing the steps outlined in Ref. 18 (Section 23), namely by 
exploiting the scaling equation featured by scaling func- 
tions and wavelets, Eqs. (9-10). Thanks to the localiza- 
tion of scaling functions, the number of non-trivial g^ 00 ''s 

amounts to just a few: c^ 00 ' ^ \ j\ < 2m — 2, where 
m is the order of the Daubechies wavelet family. The 
same holds for the other d vectors, which, as a further 
consequence of Eqs. (9-10), can be computed as follows: 



f)- 2 V/! a d m 

1 1 . v 

d {10) -2Vo h d i00) 



2i-\-fx—u ' 



(A6a) 
(A6b) 
(A6c) 



where {/i^}^-™.™], {fl^eli-m,™] arc the low " and 
high-pass filters, respectively (cf. Sec. II). 

The entries of the vector u^ 00 ) can instead be computed 
as explained below, where M = (fj, + v)/2, N = (fj, — u)/2, 
and the scaling relation Eq. (9) is used: 



it' 



(00) 



/d 
dx cj)(2x - J - (i) x— <f>(2x + J - v) 

/d 
dx <j>(2x -J-M-N) x— (j>(2x + J - M + N) 
jw,jv x 

/d 
dx cp(x — J - N) (x + M)— </>(x + J + N) 



M,N 



M,N 



,(00) 



,(00) 



W 2(J+N) + M< ^2(J+N) 



(A7a) 
(A7b) 
(A7c) 
(A7d) 
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Integrating by parts Eq. (A4), it is easy to prove that 

wf 0) ' (11) = -/f 1 '-^; (A8a) 



(01) (10) 

w j = — -j > 



(A8b) 



and, as a corollary, that Wq = —1/2 = Wq . Moreover, 



\j\ > 2m — 2 =>■ w 



(00), (01), (10), (11) 
3 



0. Restoring the 



original indexes fj,, v in (A7d), and rearranging terms in 
a convenient way, we are left with 



(A9) 



where 



and 



Aj,k = Sjm - '^h^hv Sp-v+zj,* (A10) 



h^^^CUj- (All) 



Namely, the entries of u^ 00 ' can be found as the solution 
of the linear system of equations (A9), which is neither 
under- nor over-determined. Formally equivalent linear 
systems, although with different numerical coefficients, 
lead to the solution for the vectors w^ 01 \ w^- 10 \ w^ 11 ' . 
For instance, i// 01 ) is the solution of 



with 



and 



i ,„( 01 ) - u 



w 



(A12) 



A 'j,k = fi.J.k ~22 h i x9v <W^+2J,fc ( A 13) 
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